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Abstract Reptiles are the most morphologically and 
physiologically diverse tetrapods, with the squamates 
having the most diverse habitats. Lizard is an important 
model system for understanding the role of visual ecology, 
phylogeny and behavior on the structure of visual systems. 
In this study, we compared three opsin genes (RH2, LWS 
and SWSI) among 49 reptile species to detect positively 
selected genes as well as amino acid sites. Our results 
indicated that visual opsin genes have undergone divergent 
selection pressures in all lizards and RH2 and LWS suffered 
stronger positive selection than SWSI. Twelve positively 
selected sites were picked out for RH2 and LWS. Moreover, 
many diagnostic sites were found between geckos and 
non-gecko lizards, most of which were located near the 
positively selected sites and some of them have already 
been reported to be responsible for significant shifts of the 
wavelength of maximum absorption (Amad. The results 
indicated that the gecko lineage accelerated the evolution 
of these genes to adapt to the dim-light environment or 
nocturnality as well as the switch between nocturnality and 
diurnality. 
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1. Introduction 


Reptiles are the most morphologically and physiologically 
diverse tetrapods, which diverged from early tetrapods in the 
late Carboniferous period approximately 310-320 million years 
ago along with adaptive evolution (Hughes, 2008). In particular, 
squamates represent the most diverse radiation of terrestrial 
vertebrates and occupy ecologically diverse habitat niches. Even 
though, most of them rely heavily on vision systems as their 
primary sensory modality which is closely related to the ability 
and means to interact with their environments (Fleishman et al., 
2011). Vision in vertebrates is especially well studied, and studies 
of the evolution of their visual pigments have been able to both 
identify evolutionary changes, and to ascribe such changes to 
adaptive evolutionary processes (Hughes, 2008). 

The vertebrate retina contains two types of photoreceptor: 
rods which are responsible for dim light and mediate scotopic 
vision, and cones acting in daylight, which mediate photopic 
vision and are necessary with color vision (Bowmaker and 
Hunt, 2006). Comparative studies across all of the major 
vertebrates groups have detailed that there are one rod pigment 
RH1 (rhodopsin) and four cone pigments encoded by distinct 
opsin gene families: LWS/MWS (long/middle-wavelength 
sensitive opsin gene) maximally sensitive in the red-green 
spectral region from 490-570 nm, RH2 (RH 1-like) sensitive in 
green from 480-535 nm, SWS2 (short-wavelength sensitive 2 
opsin gene) sensitive in the blue from 410-490 nm and SWSI 
(short-wavelength sensitive 1 opsin gene) sensitive in the violet- 
ultraviolet from 355-440 nm (Yokoyama et al., 2000). 

However, different lifestyles were in harmony with 
different habitat circumstance and various kinds of visual 
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system. For instance, molecular evidence implies that LWS was 
inactivated on at least five cetaceans which are known to adapt 
to dim-light field dominated by blue light (Meredith et al, 2013). 
Several groups including cetaceans, pinnipeds, bats, and mole 
rats, inhabit in low photopic environments, some of which 
have lost their SWSI opsins (Zhao et al., 2009). Furthermore, 
the significant divergent selection in bats SWSI and LWS was 
reported to be occurred in response to echolocation ability and 
foraging habitat, respectively (Gutierrez et al., 2018). Among 
Squamata, dedicated fossorial snakes lost all visual opsins other 
than RHI, whereas all other snakes (including less dedicated 
burrowers) also have functional SWS1 and LWS opsins. In 
contrast, all lizards, even highly fossorial amphisbaenians 
express functional LWS, SWS1, SWS2, and RH1 genes, and most 
also express RH2 (Simões et al., 2015). 

Lizards provide an important model system for 
understanding the role of visual ecology, phylogeny and 
behavior on the structure of visual systems in Squamata 
(Fleishman and Persons, 2001; Loew et al., 2002; Stuart-Fox et 
al., 2007; Macedonia et al, 2009). This extensive group includes 
true chameleons, iguanid lizards, monitor lizards, geckos, and 
skinks. Some of them are diurnal, and some crepuscular or 
nocturnal members. Possible changes during transformation 
of the lifestyle in these reptiles from diurnal to nocturnal have 
induced some adaptive evolution (Bowmaker, 2008). Therefore, 
it's notable to explore the evolution of opsin genes owing to 
lifestyle transformation and understand molecular features and 
specializations in squamates. For example, Anolis carolinensis, 
a diurnal lizard representive, possesses a complete set of opsin 
paralogs, however, a nocturnal gecko, Gekko japonicus only has 
three functional opsins: SWSI, LWS, and RH2 (Liu et al., 2015). 
This result from genomic sequence data is consistent with the 
inference of gecko opsin composition based on the wavelength 
) values (Kojima et al., 1992, 
Loew, 1994) and the experiments which failed to detect any 
hybridizing band using the chameleon RH1 and SWS2 cDNA 
clones as probes in the gecko genome (Shi et al., 2001). With the 
loss of RHI, the question of which functional opsin is responsible 


of maximum absorption (A 


max. 


for nocturnal vision in geckos needs to be solved. The amino 
acid replacement F89C of RH2 was suggested to allow geckos to 
receive more light (Liu et al, 2015). However, the evolutionary 
pattern of opsin they proposed was based on genome data 
of few species, and more genome data or molecular data 
from other related species are necessary to determine how 
particular changes occurred in these opsin genes. In this study, 
we aimed to further investigate adaptive evolution of gecko 
opsin genes during their behavior shift between different photic 
environments by comparing opsin genes among diverse lizard 


reptiles. 


Yao CAI et al. 
Molecular Evolution of Visual 
Opsin Genes in Geckos 


281 


2. Material and Methods 


2.1. Taxon Sampling Samples were collected by fieldwork or 
obtained by academic exchanges. Taxon sampling for surveys 
of opsin genes comprised 49 reptile species including 40 lizards 
(ingroups), 6 snakes and 3 turtles (outgroups), 27 of which were 
newly sequenced and others were downloaded from GenBank 
in this study (supplementary Table SI). Sample species includes 
seventeen geckos, six scincoid lizards; eleven Lacertoidea species; 
one Anguimorpha lizard; five iguanas; and six snakes. Samples 
used in the molecular analyses are detailed in supplementary 
Table SI. 

PCR and Sequence alignment 

Genomic DNA was extracted from muscle tissues of each 
individual following standard phenol/chloroform extraction 
methods. We attempted to amplify each exon of RH2, LWS 
and SWS] visual opsin genes using primers newly designed 
according to the sequences from Gekko gecko, Phelsuma 
madagascariensis, and Anolis carolinensis (supplementary Table 
S2-S4). All fragments were amplified in 50 uL PCR reactions: 1 
pL of each primer, 10 pmol/L; 18 uL DNA template, ca. 200-300 
ng; 5 uL 10 x PCR buffer (Mg); 4 uL dNTPs, 25 mmol/L; 05 
pL Takara Taq (Takara Biomedical, Dalian, China) and 36.7 
uL ddH;O. The PCR condition was as follows: 3 min at 94°C 
followed by 35 cycles of 30 s at 94 °C, 30 s at 52-60 °C, 90 sat 72 
°C, and a subsequent 10 min final extension step at 72 °C. PCR 
products were checked by electrophoresis on 1.5% agarose gels 
and purified with AxyPrep PCR Cleanup Kit (AxyGen MA, 
USA). Sequencing was carried out bidirectionally on an ABI 
PRISM 3730 DNA Sequence. 

The newly generated sequences were examined by 
comparison with the published ones by BLAST on NCBI, and 
then assembled by Lasergene v7.1 (DNASTAR Inc, Madison, 
WI, USA). After manual editing, all sequences were aligned by 
MEGA70 (Kumar et al, 2016). Totally 431 fragments with high- 
quality were picked and concatenated to 80 opsin genes, and 
then analyzed together. The newly sequences were deposited 
in GenBank under accession numbers KU645203-KU645282 
(supplementary Table SI). 


2.2. Molecular evolution analysis The CODEML program 
from the PAML 49 package (Yang, 2007) was used to estimate 
the non-synonymous (dy) and synonymous (ds) substitution 
rates and their respective ratio (d,/d, or w) for the RH2, LWS 
and SWSI genes among squamates. The w indicates changes in 
selective pressures, where o = 1, w < 1 and w > 1 corresponds to 
neutral evolution, purifying and positive selection, respectively. 
A well-accepted squamate phylogenetic tree was used as a 
framework tree in our analysis for each gene (Rosler et al, 2011, 
Pyron et al., 2013, Zheng and Wiens, 2016). Different datasets, 
dataset IA (all lizards), dataset IB (all squa mates), dataset II (only 
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geckos), and dataset III (non-gecko lizards) were used for PAML 
analyses to minimize the effects of sampling. 

To examine the probabilities of sites under positive selection 
in three opsin genes, we used site models M8a (B and @, = 1) 
versus M8 (Swanson et al, 2003; Wong et al., 2004) implemented 
in the CODEML program of PAML 4.9 (Yang, 2007). The 
nested models were compared using a likelihood ratio test 
(LRT) with a 7’ distribution. Positively selected sites in the M8 
were identified using a Bayes Empirical Bayes (BEB) analysis 
(Yang et al, 2005) with posterior probabilities > 0.80. Positively 
selected sites were also detected by fixed effects likelihood 
(FEL) performed in HYPHY (Pond and Frost, 2005) (via the 
www.datamonkey.org web server), with the default settings 
of significance levels of 0.2. We then performed selective 
pressure detection using TreeSAAP v.32 (Woolley et al, 2003), 
which detected selection based on 31 physiochemical amino 
acid properties. All magnitude category 6-8 changes with P < 
0.05 were used as an index for the degree of radical amino acid 
substitutions and positive selection. 

To evaluate whether positive selection was restricted to 
specific squamate lineages, we used branch models (including 
free-ratio model and two-ratio model) (Yang, 1998; Yang and 
Nielsen, 1998) and branch-site model implemented in CODEML. 
The free-ratio model (M1) that assumes an independent w ratio 
for each branch was compared to the null one-ratio model 
(MO) with the same o for all branches (Yang and Nielsen, 1998). 
Two-ratio model and branch-site model require the foreground 
branches (lineages tested to be under positive selection) and 
background branches (rest of the lineages) to be defined a priori. 
The o ratios were estimated for two branch categories (geckos 
as foreground branches versus non-gecko lizards as background 
branches) in two-ratio model, and then also compared with 
the null MO model. The modified branch-site model A with w 
varying among sites and among lineages (Zhang et al., 2005) 
was tested against the null hypothesis of no selection in any 
of the foreground or background branches. Each squamate 
lineage across the phylogeny was used as the foreground 
branch, respectively, where the remaining branches were 
treated as background branches for each gene. A false discovery 
rate (FDR) correction for multiple tests was applied to the LRT 
P values for branch-site model analysis (Anisimova and Yang, 
2007). 


2.3. Sliding Window Analyses To explore the further 
heterogeneous selection pressure at RH2, LWS and SWSI opsin 
genes across squamate phylogeny, we constructed a sliding 
window analysis using the program SWAAP (Pride, 2000) 
with window size at 30 bp and step size at 3 bp, which were 
repeated for several groups of interest. Sliding window analysis 
was performed to show variation in w value (d,/d,) along each 
opsin genes in geckos, non-gecko lizards and all lizards. 
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2.4. Mapping of positively selected sites onto opsin 3D 
structures To gain insights into the functional significance 
of these putatively selected sites, we mapped these positive 
sites onto protein structures. The 3D structures of opsin genes 
under positive selection were predicted by using the homology 
modeling software provided on the I-T'ASSER server (Zhang, 
2008). The protein sequences of positively selected genes were 
derived from the common Tokay gecko (Gekko gecko), which 
were obtained from the GenBank database. 


3. Results 


3.1. Visual pigment complement We assembled a dataset 
of opsin gene sequences for 49 squamates, comprising 17 
Gekkota geckos, 6 Scincoidea lizards, 11 Lacertoidea lizards, 1 
Anguimorpha lizard, 5 Iguania lizards, 6 Alethinophidia snakes 
and 3 turtles (outgroups). We have successfully amplified three 
visual opsin genes (RH 2, LWS, SWS1) from genomic DNA of 
the lizards from 20 families. All sequences were assembled and 
concatenated correctly for the protein coding region of RH2, 
LWS, and SWSI, which contain 1065 bp (354 aa), 1110 bp (369 
aa), 1044 bp (346 aa), respectively. The sequence alignments had 
no premature stop codons or frame-shift mutations, which 
indicated the presence of functional proteins in squamates. 


32. Selection on opsin genes We found M8 that incorporated 
selection fitted the data better than the neutral model, M8a, 
at two opsin genes (RH2 and LWS, P « 0.05), suggesting that 
these genes were subjected to positive selection in squamates. 
However, for SWS], the LRT showed no significant difference 
between the models M8 and M8a implying no positive 
selection. Using M8, the most stringent model carried out 
in PAML, a small proportion of codons (0.84%-10.17%) were 
estimated to be under selection with average w values of 1292- 
2.462 at the two positively selected genes in squamates (Table 
1). 3 (0.8496), 36 (10.17%), and 3 (0.84%) positively selected sites 
were identified respectively by the BEB approach as having 
posterior probabilities > 0.80 at RH2 in dataset IA (all lizards), 
dataset II (only gecko), and dataset III (non-gecko lizards). 
Meanwhile, 6 (1.6396), 7 (1.8996), and 8 (2.17%) positively selected 
sites were identified respectively at LWS in dataset IA, dataset 
IB (all squamates), and dataset II. However, no positively 
selected sites were detected at LWS in dataset II. When we used 
a significance threshold of 0.95 for posterior probabilities, the 
number of positively selected amino acids decreased to two, 
thirteen, and two for dataset IA, II, and III of RH2, while one, 
one, and two for dataset IA, IB, and II of LWS. FEL, performed 
in HYPHY (Pond and Frost, 2005), was also used to test for 
selection in the three opsin genes. HYPHY can improve the 
estimation of the ¿ value by incorporating variation in d; 
whereas d; is fixed across sequences for all the PAML-based 
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Table1 Positively selected sited detected using two maximum likelihood (ML) methods across squamate phylogeny. 
Test of selection Sites under selection identified by ML methods 
-Ln(M8a) -Ln(M8) -2AL w value PAML M8 FEL? No. of sites’ % of sites 


IA: all lizards 8635191 8629685 11012 2.438 169, 209, 213 9 2.54% 
ka 13, 19, 24, 39, 41, 49, 50, 86, 88, 97, 100, 104, 108, 112, 24, 39, 112, 154, 
(354 aa) II: only gekkonids 3850724 3842002 17444 2111 154, 158, 166, 168, 205, 209, 210, 213, 214, 216, 259, 158, 162, 169, 

R 263, 270, 271, 273, 279, 282, 286, 290, 292, 299, 336 205, 209, 213 

Hl: non-gekko lizards 5303745 5299837 7816 2462 165,169, 213 

IA: all lizards 8911865 8909803 4124 1.307 72, 112, 175, 216, 279, 290 19, 72, 116, 178, 3 0.86% 
LWS  IE:all squamates 10231735 10228664 6.142 | 1292 31, 112, 175, 216, 279, 282, 290 186, 233, 241, 
(369 aa) l: only gekkonids ^ 4165037 4160029 10016 ^ 264 31,112, 175, 182, 188, 229, 279, 282 279, 290, 292, 

Il: non-gekko lizards 5445293 5444968 0.65(NS) 1585  . 351 E _ 

IA: all lizards 7691336 7690724 122(NS) 1447 _ 37, 51, 83, 157, ES X 
SWS1 IB: all squamates 9505075 9505.077 0.04(NS) 1916 = 213, 293 = = 
(346 aa) IE only gekkonids 3324288 3324288 O(NS) 1 E e s E 


Il: non-gekko lizards 4752117 4751884 047(NS) 1454 _ 


The positively selected sites picked out by two methods are shown in bold; ` indicates P < 005 whereas NS indicates ‘Not significant’; * Codons were detected 
by M8 model in PAML using a Bayes (BEB) analysis with posterior probabilities > 0.80; ` Codons were detected by FEL implemented in HYPHY with 
significance levels of 02; * No. of sites indicated positively selected sites identified by both ML methods. 


analysis. FEL analysis showed that significant signs of positive 
selection were detected at the three opsin genes with many 
more positively selection sites than that identified by M8 (Table 
1). Combining the two different maximum likelihood (ML) 
methods, a total of 12 positively selected sites (9 at RH2 and 3 
at LWS) were picked out (sites in bold shown in Table 1). Sites 
that were identified to be under positive selection by two ML 
methods were regarded as robust candidates for sites under 
selection. We further employed a complementary protein- 
level approach implemented in TreeSAAP (Woolley et al, 2003) 
to evaluated destabilizing radical changes at each robust site. 
The result showed that 10 of 12 sites (83.88%) have at least three 
radical changes in properties (Table S5, Table S6) 

To evaluate whether positive selection is only limited 
to particular lineages at the three opsin genes, we first used 
branch models (including free-ratio and two-ratio models) 
that allow the w ratio to vary among branches across the 
phylogeny. The LRT tests showed that relaxation of selective 
pressure was detected at RH2 (o... | = 0.296, o 
0.156) and LWS genes (o 
= 0343, o 


as foreground branch and the remaining non-gecko lizards 


geckos non-gecko lizards 7 


= 0398, o = 0.141; w 
= 0.149) when lineage of geckos was treated 


geckos non-gecko lizards geckos 


other squamates 


or other squamates were treated as background branches. 
On the contrary, the SWSI (w =0.044) was implied to be under 
very strong purifying selection (data not shown). Using the 
stringent branch-site model, we found that w was greater than 
1 at RH2 in four terminal branches leading to Cnemaspis limi, 
Phelsuma laticauda, Gekko swinhonis, and Fylinia sp. respectively, 
among which three were nested in the Gekkota clade (Figure 
1, Table S5). As well at LWS, evidences of positive selection 
were examined within lineage Gekkota. The w values were 
greater than 1 in the last common ancestor (LCA) of Gekkota, 
terminal branches leading to Gonatodes annularis, and leading 
to Hemidactylus bowringii (Figure 1, Table S5). In contrast, 


in branch-site model significant signs of positive selection 
were also examined at SWS1 in some ancestral branches 
(LCA of Gekkota, LCA of scincids, LCA of lacertids, LCA of 
Hemidactylus, LCA of Pachydactylus vansoni + Phelsuma) and 
in one terminal branch leading to Phelsuma madagascariensis 
(Figure 1, Table S5). Furthermore, each identified branch was 
detected some codons under positive selection. The results 
of branch-site models suggested more widespread positive 
selection across opsin genes in geckos than other squamates. 
This apparent difference can also be observed in the result of 
our sliding window analyses, with many peaks significantly 
greater than 1 for geckos, while the curves were uniformly 
below 1 for the negative control clades with only one exception 


(Figure 2). 


3.3. Spatial distribution of the positively selected residues 
on opsin 3D structures A total of 12 radical amino acid 
changes subjected to positive selection identified by two ML 
methods were mapped onto the 3D structures of two opsin 
genes (nine at RH2 and 3 at LWS). Opsins are G-protein 
coupled receptors (GPCRs), which contains seven a-helical 
transmembrane domains, an extracellular domain, and an 
intracellular domain to coupled (Terakita, 2005). We found that 
11 of the 12 (91.67%) positively selected sites were localized in the 


a-helical transmembrane domains (supplementary Figure SI). 


3.4. Diagnostic amino acid sites between geckos and non- 
gecko lizards We also compared the amino acid sequences of 
the three opsin genes between geckos and non-gecko lizards to 
obtain the information of diagnostic sites. As shown in Table 
2, when putting nocturnal and diurnal geckos together, 16 
amino acid residues of RH2 differed from those of non-gecko 
lizards (at residue 40, 41, 42, 43, 49, 50, 54, 56, 60, 64, 83, 85, 98, 
285, 289, and 300, on a gray background in Table 2). However, 


residues 89 and 101 in diurnal geckos represented both status 
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Table 2 Comparison of the diagnostic amino acid sites between the opsin genes of geckos and non-gecko lizards. 


Codon position 


40 41 42 43 49 50 54 56 60 64 83 85 98 285 289 300 

nocturnal L S/A FIL M V/I/A V/A L G PERO N L/V S Ww P WAY 
ond diurnal L S/T FIL M Wil, Wh, G EO N l Š Yoro y 
Non-gecko lizards diurnal V/I C V/C I S/F T I I/L/F L K D F A A A/S L 

Codon position | E Codon position 

= 220 225 281° 298 321 EEUU 11 19 56 57 76 92^ 128 149 

Gd nocturna C LF N Y nae nocturnal A F A A A/V S i Sm 
diurna © JL IP H/Y Ww diurnal  A/ST F A A A $ J Sg 

Non-gecko lizards diurnal V/I V/I Y N Y Non-gecko lizards diurnal E W V T T A F L/M 


“indicates positively selected sites by M8 model in PAML. Underlined sites are reported spectral tuning sites; * corresponding to site 277 in Yokoyama and 
Radlwimmer, 1998; ° corresponding to sites 90 and 97 in Yokoyama et al, 2006. 


Outgroup 
2:61.79 
26.80 
Gekkota 
Sphenomorphus indicus 
Scincella modesta 
Plestiodon chinensis Scincoidea 
Plestiodon elegans 
miss Melanoseps occidentalis 
Feylinia sp. 
Lacertoidea 
Anguimorpha 
Iguania 
RH2 — Serpentes 
LWS — Diurnal species & 
SWS] — Amino acid change | 


Figure 1 Evidences of positive selection across Squamata phylogeny identified by branch-site model are shown. Three genes are marked with different 
colors, ie, RH2 (red), LWS (yellow), and SWS1 (green). The w values greater than 1 for individual lineage according to free-ratio and two-ratio are shown. 
Amino acid changes identified using TreeSAAP across colored branch were indicated by blue vertical lines. The circles were on behalf of diurnal genera 
in geckos. All detailed results were listed in supplementary Table S5. 
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found in nocturnal geckos (89C, 101A) and non-gecko lizards 
(89F, 101G), indicating an intermediate status (shown on a dark 
blue background). Among all these amino acid residues, three 
sites (41, 49, and 50, with a star in Table 2) were identified as 
positively selected sites by M8 model in PAML, and most of 
the rest are surrounding the positively selected sites. Moreover, 
certain changes at three sites (49, 83, and 89, underlined in 
Table 2) have already been reported to cause significant Amax- 
shifts in various visual pigments during vertebrate evolution 
(Yokoyama et al., 2007; Takenaka and Yokoyama, 2007). 
Furthermore, the amino acid in codon position 89 has been 
reported to be an important functional determinant of rod 
and cone visual pigments (Yokoyama and Blow, 2001), and 
the amino acid replacements F89C enables RH2 to have rod 
pigment-specific biochemical characteristics (Yokoyama and 
Blow, 2001; Yokoyama and Tada, 2010). F89C was found in 
all geckos included in this study except one diurnal species 
(Gonatodes annularis). For LWS and SWSI, three and eight 
diagnostic sites can distinguish geckos and non-gecko lizards, 
respectively, and one to two sites have been reported to be 
important spectral tuning sites (Yokoyama and Radlwimmer, 
1998; Yokoyama et al., 2006). Among these sites at SWSI, in 
diurnal geckos both 85A (same to nocturnal gecko) and 85S 
(same to diurnal lizards) were detected. Besides, two sites which 


2.0 
1.5 


1.0 
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are shared between nocturnal gecko and diurnal lizards (298N 
and 321Y) were totally different in diurnal geckos (298H/Y and 
321W, respectively) (on a light blue background in Table 2). 


4, Discussion 


Most vertebrates use a combination of rod and cone 
photoreceptors to enable vision in conditions ranging from 
starlight to direct sunlight. However, geckos are known to 
possess a mix of all-rod (nocturnal geckos) or all-cone retinas 
(diurnal geckos), which are all suggested to be evolutionary 
derived from cones by both morphological and molecular data 
(Roll, 2000; Roll, 2001b; Zhang et al., 2006). Based on comparative 
retinal and photoreceptor morphology, Walls (1942) proposed 
a scenario for the evolution of photoreceptors in geckos that 
originated with the loss of rods in a diurnal ancestor of geckos, 
resulting in an all-cone retina. From this diurnal all-cone state, 
ancestral geckos transitioned to nocturnality and evolved an 
all-rod retina through photoreceptor transmutation. Some 
extant gecko lineages subsequently reverted back to diurnality 
and re-evolved all-cone retinas through transmutation of 
the "rods" back into cones (Schott et al., 2019). Recent genome 
researches provided the evidence of a RHI deficiency at the 
genome level in nocturnal geckos (Hara et al, 2018; Liu et al, 


915 (bp) 


1065 (bp) 


> habas, ` 


765 915 (bp) 


Figure 2 Sliding Window Analyses with window size at 30 bp and step size at 3 bp to show variation in w value (dy/d,) along each opsin genes (A: RH2, 
B: LWS, C: SWS1) between only geckos (Red), non-gecko lizards (Green) and all lizards (Black). 
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2015) , and the loss of RH1 occurred earlier than that of SWS2, 
which was in agreement with the hypothesis that the ancestors 
of modern geckos were diurnal lizards without rod opsin 
(Liu et al., 2015). Moreover, the F89C amino acid substitution 
in RH2 gene shared between the Madagascar ground gecko 
and Japanese gecko suggested a possible functional switch or 
compensatory change in cone pigment of geckos that allow 
them to receive more light (Hara et al., 2018; Liu et al., 2015). 
However, the evolution pattern of opsin functional switch was 
proposed based on genome data of few species. More data from 
more related species are needed. 

As well known, most lizards are strictly diurnal like 
iguanids, chameleonids, agamids, scincids, lacertids, anguids, 
pygopodids, and varanids (Roll, 2001a). However, as the largest 
extant lizard families, most geckos are nocturnal, with some 
lineages have been reversed back to diurnal. In this study, three 
opsin genes (RH1, SWSI, and LWS) from nocturnal as well as 
diurnal geckos, and several lizard species were sampled and 
analyzed to further investigate the possibly important amino 
acid changes during evolution or adaptation. 

Site models (M8) identified 9 and 3 codons as robust sites 
under positive selection by at least two Maximum Likelihood 
methods for RH2 and LWS respectively. However, model M8 
failed to identify any positively selected sites for SWS1 opsin 
(Table 1), suggesting somewhat purifying selection on SWS1 
opsin gene. In Branch-site model, w was found to be greater 
than 1 at RH2 in four terminal branches, among which three 
were nested in the Gekkota clade. Similarly at LWS, evidences 
of positive selection were examined within lineage Gekkota 
(Figure 1). Moreover, significant signs of positive selection were 
also examined at SWSI in some ancestral branches and some 
codons under positive selection were found which were not 
supported by the site model. Both of the results of branch- 
site models and the sliding window analyses supported more 
positive selection across opsin genes in gecko lineage than in 
non-gecko lizards. 

Among geckos, there are only 15 genera being diurnal (such 
as Phelsuma, Sphaerodact ylus, and Gonatodes) that have pure- 
cone retinas as a consequence of a change of their habitat (Roll, 
2001a). Previous researches revealed difference in anatomical of 
the retina between nocturnal and diurnal geckos. They pointed 
that nocturnal geckos had photoreceptors with large rod-shaped 
outer segments and colored oil droplets in the inner segment 
whereas diurnal species had much smaller somewhat more 
cone-like outer segments and colorless oil droplets (Roll, 2000; 
2001a). However, our results provided an additional molecular 
supplement of the opsin genes between nocturnal and diurnal 
geckos. In all our samples, three species (Gonatodes annularis, 
Phelsuma madagascariensis, and Phelsuma laticatus) on behalf of 
diurnal genera had undergone stronger adaptive selection than 
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others geckos resulting in transformation of diurnal ambient 
light (Figure 1). Previous researches have suggested that RH2 
opsin amino acid site 89 is an important functional determinant 
of rod and cone visual pigments (Shi et al., 2001), and the amino 
acid replacements F89C in geckos enables RH2 to have rod 
pigment-specific biochemical characteristics to allow them to 
receive more light in dim-light environment (Kojima et al, 1992; 
Shi et al, 2001; Yokoyama and Tada, 2010; Liu et al, 2015; Hara 
et al., 2018). But these researches are based on the comparison of 
few opsin sequences. In our study, comparing the orthologous 
RH2 pigment among 39 species, we deduced that all geckos 
have amino acid site 89C except G. annularis (diurnal) and all 
non-gecko lizards were 89F (Table 2). 89 F in G. annularis is 
likely to be the result of reversion due to the re-adaptation to 
daylight. Similar change was also found at site 85 in SWSI, 
which was also a reported spectral tuning site. Furthermore, we 
also found sites like 298 and 321 in LWS, where nocturnal gecko 
and diurnal lizards share the same residue; however diurnal 
geckos have different residues. These changes are likely to occur 
during the evolution of geckos from nocturnality to diurnality. 

In conclusion, our results indicate that opsin genes are under 
different selection pressures in all lizards, and gecko lineage 
has accelerated the evolution of these genes to adapt to the 
dim-light environment or nocturnality as well as the switch 
between nocturnality and diurnality. We compared amino 
acid sequences of three opsin genes (RH2, SWSI, and LWS) 
from representatives of many extant lizard families, which 
provide clues and evidences for further study of the molecular 
specialization and evolution to response to different lifestyle 
(nocturnal or diurnal) in geckos. 
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Appendix 


Figure S1 The predicted three dimensional structures (ribbon diagram) of RH2, LWS, and SWSI. a-helices, Q-sheets, and loops are represented by spirals, 
arrows, and lines, and are coloured in different colours, yellow, and red, respectively. A: RH2; B: LWS; C: SWSI. The positive selection sites (24, 39, 112, 
154, 158, 169, 205, 209, 213 for RH2; and 72, 279, 290 for LWS) are marked with red globes. 


Table S1 Samples with taxon and Genbank accession number. Sequences newly generated for this study are indicated in bold. 


Clade/Higher : Accession Numbers 
Taxon Spone RH2 LWS SWSI 
Hemiphyllodactylus yunnanensis KU645203 KU645230 KU645257 
Hemidactylus bowringii KU645204 K-U645231 KU645258 
Hemidactylus frenatus KU645205 KU645232 KU645259 
Gekko auriverrucosus KU645206 KU645233 KU645260 
Gekko hokouensis KU645207 KU645234 KU645261 
Gekko scabridus KU645208 KU645235 KU645262 
Gekko swinhonis KU645209 KU645236 KU645263 
Gekko japonicus KU645210 KU645237 KU645264 
aad Gekko chinensis KU645211 KU645238 KU645265 
Gekko gecko M92035 M92036 AY024356 
Phelsuma madagascariensis AF074044 AF074043 AF074045 
Phelsuma laticauda KU645212 KU645239 KU645266 
Cnemaspis limi KU645213 KU645240 KU645267 
Ptenopus garrulous KU645214 KU645241 KU645268 
Pachydactylus vansonic KU645215 KU645242 KU645269 
Phyllodactylus victus KU645216 KU645243 KU645270 
Gonatodes annularis KU645217 KU645244 KU645271 
Plestiodon chinensis KU645218 KU645245 KU645272 
Plestiodon elegans KU645219 KU645246 KU645266 
Scincoidea Melanoseps occidentalis KR336743 KR336713 KR336718 
(6 species) Feylinia sp. KR336742 KR336714 KR336717 
Sphenomorphus indicus KU645220 KU645247 KU645273 
Scincella modesta KU645221 KU645248 KU645274 
Amphisbaena infraorbitale KR336730 KR336704 KR336719 
Amphisbaena alba KR336729 KR336705 KR336720 
Amphisbaena sp. KR336728 KR336706 KR336721 
. Bachia flavescens KR336731 KR336703 KR336715 
T Takydromus septentrionalis KU645222 KU645249 KU645275 
Takydromus sexlineatus KR336727 KR336707 KR336722 
Podarcis sicula AY941829 -- -- 
Eremias multiocellata KU645223 KU645250 KU645276 
Eremias velox KU645224 KU645251 KU645277 
Eremias przewalskii KU645225 KU645252 KU645278 
Eremias arguta KU645226 KU645253 KU645279 
Anguimorpha Ophiodes striatus KR336732 KR336708 KR336716 
Leiolepis reevesi KU645227 KU645254 KU645280 
. Phrynocephalus vlangalii KU645228 KU645255 KU645281 
een Calotes versicolor KU645229 KU645256 KU645282 
Anolis carolinensis NM001291394 U08131 XM. 003229162 
Uta stansburana DQ100324 DQ129869 DQ100325 
Tropidophis feicki -- KR336709 KR336723 
Xenopeltis unicolor -- FJ497235 FJ497234 
Serpentes Python regius == FJ497238 FJ497237 
(6 species) Polemon collaris -- KR336710 KR336724 
Pseustes poecilonotus s= KR336711 KR336725 
Atractus flammigerus -- KR336712 KR336726 
Pelodiscus sinensis XM_006119345 -- XM_006110204 
Testudines 


(3 species) 


Chrysemys picta bellii 
Chelonia mydas 


XM_005309675 
XM_007063507 


XM_005281282 


XM_005281289 
XM_007067421 


Table S2 Primers for RH2 exons amplification from genomic DNA. 


Exons Primers Sequence 5’-3’ 


RH22-F1 5'-GTCAAGTTGCACTCTGGTC-3' 


Exon2 
RH22-R1 5'-CTGGACCAACCAAACAG-3' 


RH24-FI 5'-GCAGCTGCTCAGCARCAGGA-3' 


Exon4 
RH24-R1 5’-CTGTTTGTTCAWGAGGACATA-3’ 


Table S3 Primers for LWS exons amplification from genomic DNA. 


Exons Primers Sequence 5’-3’ 


LWS2-F1 5’-GCCCTTTCGAWGGTCCAAACTAT-3’ 


Exon2 
LWS12-R3 5'-CGCAAGAGGATACGACGTACC-3' 


LWS4-F1 5'-GTACTGGCCYCAYGGTCTGAAGA-3' 


Exod 
TOM LWS4-R1 5-CGCACGRATMGCCAMCCACACTTG-3' 


LWS6-F1 5'-TTCCGTAACTGCATAWTGC-3' 


Bond 
ENT LWS6-R2 5'-TTATGCYGGTGATACAGAGGA-3' 


TableS4 Primers for SWSI exons amplification from genomic DNA. 


Exons Primers Sequence 5' - 3' 


SWS12-F1 5'-GBCTGGTCACYGGCTGGTC-3' 


Exon2 
SWS12-R1 5'-CTGCTCCACCCGAAAGAGGG-3' 


SWS 14-F2 5’-GTAGCCGCYCAGCARCAGGAGT-3’ 


Exon4 
SWS14-R1 5'- CTGTTTGTTCATGAAGCA-3’ 


Table S5 Positive selection detected base on branch-site model in all squamates. 


Genes lineages Models -InL 2AlnL P values o values Positively selected sites 


1661 0.533, 230S 0.731, 275V 0.899, 280A 


Ma0  10298.420 w0=0.03,01=1,02=1 


e. G.annularis M 465 0031 0.994**, 292V 0.598, 296A 0.775, 298N 0.739, 
a 10296.095 40-003,)1-1,02-8.593  331Y 0.99575 
Ma0 10307.880 0=0.032,91=1,92=1 ; I i 
DWS ¿upaqa a UE ol-lo 157F 0.817, 1661 0.687, 170S 0.810, 186 A 
Ma 10303.666 40-0.031,01-1,02-999.00 0.513 
Ma0 10305.616 0=0.031,01=1,02=1 
g. Gekkota à 4166 0041 ^" — 217E 0.654,220C 0.990**, 281F 0.670 


Ma 10303.533 4020.031,0121,92-26.797 


Note: Sites with P > 0.95 and w2 >1 are indicated in bold. Any P value higher than 0.99 was designated with two (**) asterisks, and P value 
between 0.95 and 0.99 was designated with one (*) asterisk. 


Table S6 Amino acid sites identified as positively selected sites detected using PAML, Datamonkey and TreeSAAP at RH2, LWS and SWS1 


opsin genes. 
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(Continued Table S6) 


PAML (P > 08) Datamonkey* TreeSAAP* 
AN. Branch-site SLAC FEL REL 
site model odel (P<02)  (P<02) (P >50) 


Amino Acid Properties AA change 


à 
= 
š 


PPP SoA 


Š 


N.Bu 


PPP ASS 


x 


PPP 


P.N eK°pK'a, 


i 


PBPFRP 


205 0983 0.117 0.096 482.93 pK IM 


209 


Š 


0171 0.186 16595 N.BpK'R,H,H, 


pK’ MBI 


R 
E 
L 


5 


210 0911 D,BP.FR, 


X 
t 


E 
E 
L 


> 
E 
- 
L 


l: 
+ 
= 


* 
< 
L 


z 
E 
< 
f 


* 
£ 


251 0.893 


263 0.953 


Baci 
< 
i 


27 0972 


> 
> 
5 


a 
za 
i 


^ 
"nj 
$ 


282 0919 K°P pE,P 


I 


"S 


ED 


(Continued Table S6) 


PAML (P > 08)" Datamonkey* TreeSAAP* 
Gene AA , Branch-site SLAC FEL REL p ; : 
Site model moda. (P « 02) (P « 02) (P >50) Amino Acid Properties AA change 
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(Continued Table S6) 


PAML (P > 08)" Datamonkey* TreeSAAP* 
Gene AA , Branch-site SLAC FEL REL p ; : 
Site model moda. (P « 02) (P « 02) (P >50) Amino Acid Properties AA change 
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PAML (P > 08) Datamonkey* TreeSAAP* 
Gene AA* B Branch-site SLAC FEL REL 
pened odel (P«02  (P<02) (P >50) 


Amino Acid Properties AA change 
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Note: a: Codons identified by both M8 model and FEL were in bold; b: Codons identified by the site model and branch-site model in PAML 
with P > 0.8; c: Posterior probabilities detected in Datamonkey of which SLAC, FEL and REL; d: Amino Acid Properties and AA change were 
calculated in TreeSAAP. 


